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CHAPTER I 


INTRODUCTION 

During the last few years, much attention has been di- 
rected toward thunderstorm-related approach and departure 
accidents involving civil and military jet transports. Vari- 
ous investigations of these accidents attribute the cause to 
severe wind shear. In particular. National Transportation 
Safety Board-NTSB (1976a) cited such a cause for the tragic 
end to Eastern Airlines Flight 66, at JFK Airport in New 
York, on June 24, 1975. Fujita and Byers (1977) and Fujita 
and Caracena (1977) have examined this and other accidents, 
arriving at the so-called "downburst-spearhead echo" theory 
of accident cause. 

There have been a number of other reports related to 
the aircraft wind shear problem. Pertinent additional NTSB 
reports are NTSB (1974, 1976b). Fundamental studies describ- 
ing how an aircraft responds to fluctuations in the wind are 
few. The authors of this report have produced McCarthy and 
Blick (1976), McCarthy ^ ah. (1978a, 1978b). Recently Frost 
and Crosby (1978) and Frost and Reddy (1978) also examine the 
relationship between wind shear and aircraft response. 

All too often, wind shear is examined without relating 
the shear to aircraft response, as in Fujita and Caracena 
(1977) . Basic response theory of aircraft is examined in 
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standard texts such as Roskam (1972) and McRuer et (1973), 

without relating it to real atmospheric winds. We describe a 
study where actual winds measured by aircraft, in the thun- 
derstorm environment, are used as input to a numerical simu- 
lation of aircraft flight, to determine response, we have 
posed a number of questions regarding aircraft landing or 
taking off in the thunderstorm environment: 

1. Can the real winds near the earth's surface, which 
might adversely affect flight, be obtained? 

2. Can a numerical model be developed which can test 
real wind profiles to determine whether an approach 
will be dangerously deteriorated? 

3. With an affirmative answer for questions 1 and 2, 
can we obtain a reasonable estimate of whether a 
particular aircraft will experience difficulty, 
knowing the aircraft-encountered winds? 

c. C m a real-time detection and warning system be 
developed using the information acquired in the 
pursuit described above? 

With. certain limitations, and yet a few uncertainties, we 
believe all of the above questions are answered in the posi- 
tive . 

To answer (1) , we were fortunate to obtain several hours 
of three-dimensional wind data, collected by a meteorologi- 
cally instrumented Queen Air airplane, on loan from the 
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National Center for Atmospheric Research (NCAR) , which was 
flown within 500 m of the surface in the Oklahoma thunderstorm 
environment. Additional wind data were obtained from the 
National Severe Storms Laboratory (NSSL) instrumented tall 
tower in Oklahoma City. Finally, subjective wind profiles 
were obtained for an examination of Eastern Flight 66 crash 
from NTSB (1976a) , data apparently described in Fujita and 
Caracena (1977) . 

A relatively simple numerical model for aircraft re- 
sponse was developed using a standard set of aircraft equa- 
tions of motion, to pursue question (2) . The model concen- 
trated on fixed-stick approaches, or in other words, no pilot 
control inputs were considered. 

Once we obtained the wind data and the model, we con- 
ducted approximately one thousand simulations of aircraft 
response to wind shear conditions associated with thunder- 
scorms, including a graphic simulation of Eastern Flight 66. 

It was possible to arrive at a parameter we called the Ap- 
proach Deterioration Parameter (ADP) , which can be used to 
determine quantitatively whether a particular airplane can 
safely approach in specific wind shear conditions. 

Although our work to this point has not been used as 
a real-time wind shear detection and warning system, we be- 
lieve that our technique makes such a system feasible. Plans 
are underway to test this contention. 
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CHAPTER II 


THE MODEL 

1. The Equations of Motion 

The equations used in the present analysis are essen- 
tially the classical equations of dynamic stability available 
from the literature, after Roskam (1972) and McRuer al . 
(1973) , with certain refinements made to describe the speci- 
fic system under consideration. The refinements to the clas- 
sical equations include the addition of wind gust forces. 

The aircraft trajectory model employed in this study 
was derived based on the following assvmiptions : 

a) The earth is flat and nonrotating. 

b) The acceleration of gravity, g, is constant 
(9.81 sec ^) . 

_3 

c) Air density xs constant (1.225 kg»m ). 

d) The airframe is a rigid body. 

e) The aircraft is constrained to motion in the ver- 
tical plane. 

f) The aircraft has a symmetry plane (the x-z plane) . 

g) The mass of the aircraft is constant. 

h) Once the aircraft is trimmed, its throttle setting 
and elevator deflection angle are not changed 
(fixed-stick assumption) . 

i) The aerodynamic stability derivatives are constant 
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within the altitude and Mach number range experienc- 
ed in this investigation. 

The longitudinal linear differential equations of mo- 
tion referenced to stability axes can be expressed as 

Cl - X^u - X^w + g0cos0j^ = 

- Vg - Vg 


-Z^u + w 


Z w - Z.<!? - Z 0 + g0sin0- 
w w q ^ 1 


- 


Z u 
u g 


- Z w 
w g 


( 2 ) 


0 - M u - M w 
u w 


w 


»q6 = 


«6®E - “u“g 






M w 
w g 


(3) 


The auxiliary relationship needed to convert the motion 
variable of Eqs. (l)-(3) into altitude is. 


h = -w + Uj^0 . 


(4) 


We selected to model a medium-size jet transport air- 
plane (Boeing 727 class) because of its high use frequency, 
and its involvement in the Eastern 66 thunderstorm-related 
crash. 

Eqs. (l)-(4) can be integrated numerically by Runge-Kutta 
schemes to obtain time-history, motion as the aircraft flies 
along the glide slope or they can be used to develop Bode 
plots showing the variation of altitude and velocity 
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amplitudes for input sine wave gusts (horizontal and vertical) . 

Table 1 contains a list of the airplane characteristics 
used in the simulation. 


Table 1. BOEING 727 CLASS - Medium size jet transport char- 
acteristics in landing mode. 


Mass of Airplane 63,958 kg (141,000 lb) 

Mass Moment of Inertia, I^. . . . 6. 1 x 10° kg-m'^ (4 . 5 x 10° slug-ft”^) 


wing Area 

145 

; m^ (1560 ft^) 

Fowler Flaps 6_ 

X!j 


Glide Angle, 9^ 


deg 

Initial Velocity, 

72 

m* sec ^ (140 kt) 

-1 


. , -3 -1 - 

X = -0.04065 sec 
u 

M 

w 

= -7.04 X 10 m • sec 

-1 

^ -4 -1 

X = 0.0738 sec 

M- 

= 2.69 X 10 m 

w 

w 


Z = -0.27263 sec"^ 


— -0.3228 sec ^ 

u 

-1 

q 

= 0 

Z = -0.622 sec 


w 

6 


Z. = -0.0257 

Z^ 

= -2.675 m-sec ^*rad ^ 

w 

6 


Z = -2.44 m-sec”^ 

q 


= -0.503 sec ^-rad ^ 


M 

u 

= 0 


2. Bode Plots 

In order to obtain some insight into the effects of 

horizontal and vertical wind gusts on aircraft speed and 

altitude. Bode plots can be constructed from Eqs. (l)-(4) for 

h/u , h/w , u/u and u/w , where for example, u/u represents 
9 9 9 9 9 
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the perturbation airspeed u response to a perturbation gust 
in the longitudinal component u (a u/u of 10.0 means, for 
example, that a 1 m*sec~^ longitudinal gust input gives an 
airspeed output response of 10 m- sec ^) . 

In order to obtain the Bode plots it is necessary to 
obtain the transfer functions by first taking the Laplace 
transform of Eqs. (l)-(4) and then solving for the ratios 
of the Laplace transform variable. When this is done, one 
obtains the following equations; 


u 

w 


1 

D 


-X. 

w 


-Z. 

w 


-X 




w 




-(M*E 


gcos0^ 

-(Uj^+Z )S +gsin0^ 


s - V 


(5) 


h 


w 


g 


+ 



s-x 

u 


qcos 0 

-z 

u 

-^w 

- (U^ + Zg) S +gsin0 

-^u 


-MgS 

s -X 

u 

-X 

w 

-X 

w 


-Z 

u 


-^w 


-^u 


-Mg/V3_)S +M^] 
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u 

u 


_1 

D 


h_ 

u 


1 

DS 


-^u 


■\j gcose. 

L 

-z 

u 

S(1 - 

^w^ -^w "(Vj^ +Zg)S + gsine 

-^u 


S-X 

u 

-X 

u 

gcos0^ 


-^u 

-^u 

-(Vj^ +^g)S +gsin0j^ 


-Mu 

-M 

u 

s' - M S 



( 7 ) 


DS 


s-x„ 

-X 

-X 

u 

w 

. 1 

-Z 

S(1-V 

-z 

u 

1 

-Mu 


-M 

1 


u 


u 


u 


( 8 ) 


where 


D = 


S-X. 


-Z 


u 


-M. 


u 


-X. 


w 


-z„ 


-(M*s +M„) 


geos 9 

-(V, +Z )S +gsin0. 


S -M S 

q 


(9) 


Now to obtain the Bode plots, S is replaced in Eqs.(5)-(8) 
by i(X) (where i is the imaginary number . Ihese equations 

can now be manipulated to obtain the magnitude ratios of the 
variables and the phase angle between them. These relations 
are plotted in Fig. 1. 
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Fig. 1. Wind gust Bode plots for Boeing I'll class airplane. 
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It is interesting to note that a horizontal gust pro- 
duces a large peak in the aircraft velocity perturbation u 
and a lesser peak in altitude perturbation h, at the aircraft 
angular phugoid frequency of 0.164 rad- sec The peak re- 

sponse for aircraft velocity perturbation is close to 20 db! 
This means that if a steady sinusoidal horizontal ' gust input 
of 2.1 m-sec ^ (4 kt) were encountered, the aircraft would 
respond with a sinusoidal velocity perturbation of approxi- 
mately 21 m-sec ^ (40 kt) ! This would mean that at one point 
in its cycles the aircraft would approach the stall speed of 
the aircraft, 51.5 m-sec ^ (100 kt) , and a maximum speed of 
92.7 m-sec~^ (180 kt) at another point during each sine wave 
cycle. Fig. 1 also indicates the aircraft perturbation alti- 
tude is also significantly affected by horizontal gusts when 
frequencies are at or below the phugoid frequency. Vertical 
sinusoidal gusts do not affect the aircraft velocity as much 
as the horizontal gusts. However, vertical sinusoidal gusts 
do produce larger changes in the airplane altitude then hori- 
zontal sinusoidal gusts at angular frequencies below 0.08 
rad -sec 

It is interesting to note that the height of the re- 
sonant peaks which occur at phugoid frequency in Fig. 1 can 
be related to several aircraft characteristics. When the 
characteristic equation of Eqs.(l)-(3) is determined (Eq. (9) 
set equal to zero), it can usually be written as the product 
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of two quadratic factors as sTiown by 


(s + 2C 


ph^ph 


s + 


2^ 


"'ph 


/ 2 
(s + 


2C U) s 
■sp sp 


+ U) 


sp 


') = 0 


( 10 ) 


Each of the two quadratic factors contributes to the 

total Bode plots shown in Fig. 1. Usually the long period 

(phugoid) damping factor is much smaller than the short 

period damping factor Q . It can be easily shown that the 

sp 

resonant peak on a Bode plot due to a quadratic factor is 
-20 log^Q (2^). Hence for ^ less than 1.0 the resonant peak 
has a positive value. As Q approaches zero, the peak becomes 
very large. It can be shown that for subsonic aircraft 
can usually be approximated by the following: 


^ph 


_1 ^ 
■{2 


( 11 ) 


It can also be shown that the exponential phugoid time 
constant, T, is inversely proportional to Gph^ph* Since 0 )^^^ 
can be approximated by 1.4g/u^, then the phugoid mode time 
is approximated by 


T = 


^1 


(12) 


Hence Eqs . (11) and (12) indicate that low values of 
and conversely large time constants which cause large 
excursions in aircraft perturbation speed (u) and altitude 
(h) are proportional to the lift-to-drag ratio of an air- 
plane. Jet transports generally have large values of 
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lift-to-drag, even witli tlieir gear down in the landing mode. 

In addition, jet transports land at a high speed (compared to 

general aviation aircraft); hence the product of U, and C^/C 

X Xj o 

combine to produce a large value of the time constant in 
Eq. (12). 

These factors would seem to suggest that light aircraft 
would not experience the same difficulty as a medium size jet 
transport when attempting to land in winds with large hori- 
zontal spectral components in the phugoid frequency range. 

it is interesting to note that three minutes before 
Eastern Flight 66 (Boeing 727) crashed at New York's Kennedy 
Airport on June 24, 1975, a light aircraft (Beechcraft Baron) 
made a successful landing although it did experience a heavy 
sink rate and an airspeed drop of 10.3 m*sec~^ (20 kt) , (from 
Fujita and Caracena, 1977) . Reconstructions of the horizontal 
wind gusts from NTSB (1976a) showed the presence of a hori- 
zontal gust that was approximately one-half sine wave in 
length of angular frequency 0.16 rad-sec”^ (close to phugoid 
frequency of a Boeing 727) at an amplitude of approximately 
13 m? sec This case will be addressed in Chapter III. 

That one. does not need an endless string of horizontal 
wind gusts at or near the phugoid frequency to cause landing 
problems can be seen in Fig. 2a. The equations of motion 
(l)-(4) were solved for a Boeing 727 class airplane (Table 1) 
entering a -3 deg landing glide slope at 500 m altitude. The 
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IAS~MSEC"'(KT) 



GROUND DIST.-KM GROUND DIST.-KM 

pig. 2. Response of Boeing 727 class airplane to (a) half sine wave 
tailwind gust and (b) elevator step input during landing. 




equations in stick-fixed mode were solved by the Continuous 
System Modeling Program (CSMP) method on the University of 
Oklahoma digital computer. A half sine wave tailwind of 
10 m- sec ^ amplitude at the phugoid frequency was assumed to 
act on the airplane at an initial altitude of 500 m. The 
half sine horizontal tailwind lasted about 19 seconds. The 
resulting aircraft velocity deviations from an indicated air- 
speed of 72.1 m-sec ^ (140 kt) were as large as 13.9 m-sec~^ 
(27 kt) . Even 76 sec after the gust had diminished to zero an 
airspeed deviation as large as 7.7 m*sec ^ (15 kt) was com- 
puted! The altitude deviations from the glide slope were as 
large as 100 m, 80 sec after the gust had diminished! The 
aircraft touched down 600 m short of the runway. So even a 
half sine wave gust near the phugoid frequency is extremely 
dangerous . 

The problem for the pilot is how to control this roller 
coaster effect seen in Fig. 2. Suppose the pilot noticed the 
decrease in indicated air speed during the first five seconds 
(Fig. 2a) and decided to correct for this by a sudden deflec- 
tion of the elevator 0.1 radian down, in order to pitch the 
nose down and pick up some speed. Fig. 2b shows that the 
airplane speed would not respond fully [IAS = 100.4 m-sec ^ 
(195 kt) ] until a lag time of about 20 sec, but by now the 
airplane response to the one-half sine wave horizontal gust 
is 86.0 m-sec ^ (167 kt) . Adding the two perturbation speeds 
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of (100.4-72.1) plus (86.0-72.1) one obtains an aircraft 
speed of 72.1 + 28.3 + 13.9 = 114.3 m- sec“^ (222 kt) . Hence 
the pilot, in making what appeared to be an obvious correc- 
tion to the elevator angle, actually made the system worse, 
i.e., he is making the aircraft unstable! Situations similar 
to the situation described in this paragraph may very well 
have occurred to aircraft in the past and caused the aircraft 
to crash during landing. It should be noted that one other 
obvious control the pilot has at his disposal is his throttle. 
However, this too could cause similar problems due to the 
long lag times required to "spool-up" the jet engines. 

Bode plots (see Fig. 3) were also plotted for altitude 
response to elevator angle and aircraft perturbation velocity 
to elevator angle. Again peaks in amplitude plots are seen 
near the phugoid frequencies. Notice, at frequencies lower 
than the phugoid frequency, the phase angle for h/s^ is about 
-180 deg (which is desirable since down elevator results in 
decreasing height) . However, just above the phugoid fre- 
quency, h and 6^, are in phase which is undesirable (down 
elevator causes the airplane to climb). The Bode plot 

indicates that u and are in phase (desirable) at frequen- 
cies below phugoid frequency but are about -180 deg out of 
phase (undesirable) above the phugoid frequency. Fig. 3 
then illustrates the problem a pilot has in adequately con- 
trolling a Boeing 727 class aircraft with the elevator in 
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■ 0 > .1 1.0 

(j rad sec~' 

Fig. 3. Elevator angle Bode plots for Boeing 727 class air- 
plane. 
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turbulent winds with strong spectral components in the phugoid 
frequency range . 

3. Spectra of Wind Fields 

A simple approximation to real wind fields is a set of 
sine waves. This is not realistic in terms of what happens 
in the atmosphere but by using such "pure tones" in the air- 
craft approach simulation model, discussed earlier, prelimi- 
nary conclusions about aircraft performance can be obtained. 
For simulated Boeing 727 approaches from 500 m along a -3 deg 
glide slope, pure sine wave inputs, at various frequencies, 
have been used. Fig. 4 shows the root-mean-square values of 
airspeed and height deviations from the normal airspeed of 
72.1 m-sec ^ (140 kt) and glide path due to a continuous 1.0 
m*sec amplitude input sine wave oscillation. Inputs in 
both the longitudinal (^g) and vertical (w^) components have 
been tested. It is obvious from the figure that maximum air- 
craft airspeed deviations occur for longitudinal inputs at 
the phugoid frequency. Vertical sinusoidal inputs do not 
affect the aircraft airspeed as much as the horizontal gusts. 
However, vertical inputs can cause much larger values of ^^h ' 
than longitudinal inputs and the maximum value occurs at about 
one-sixth of the phugoid frequency. 

since pure sine waves seldom exist in the atmosphere, 
a better way to characterize the input wind fields is in terms 
of a variance or energy spectrum. This is given by energy 
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Fig. 4. RMS altitude and velocity deviations of a Boeing 
727 class airplane to continuous one m-sec“l horizontal or 
vertical sine wave winds during a landing from 500 m alti- 
tude to the ground. 
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density as a function of frequency. Since energy densities 
are used, a frequency interval is involved when actual 
energies (or associated amplitudes) are discussed. 

There are several methods of obtaining variance spectra. 
We have chosen a technique described by Ulrych and Bishop 
(1975) as applied by Goerss and Koscielny (1977). This tech- 
nique, referred to as the "maximum entropy" or auto-regressive 
modeling method, fits an optimum auto-regressive process to 
the given data to obtain energy density as a function of fre- 
quency. In all cases, one second data are used so that the 
frequency range is 0.0 to 0.5 Hz. We have used a frequency 
resolution of 0.002 Hz in finding the energy density, so the 
energy density is calculated for 251 frequencies. 

Energy density spectra have been calculated for both 
the longitudinal and vertical wind components. For the longi- 
tudinal wind inputs, additional spectra have been produced. 
These spectra include not only the longitudinal wind input 
energy density, but also include an estimate of actual air- 
craft response; in particular, the airplane transfer func- 
tion (u/u ) is used. The energy density spectrum is multi- 
9" 

plied by the transfer function squared to give a new energy 
spectrum, which we have called the Aircraft Response to Wind 
Spectrum (ARWS) . The ARWS along with the longitudinal and 
vertical wind input spectra will be used in the discussions 
that follow. 
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CHAPTER III 


AIRCRAFT APPROACH SIMULATIONS 

1. Theory 

At the heart of our research has been the modeling of 
aircraft performance to simulate aircraft approaches in vari- 
ous input wind fields. The basis of the model and the appro- 
priate equations have been discussed in the previous chapter. 
The aircraft has been simulated to fly in a fixed-stick mode, 
in that no control or thrust changes are made during a simu- 
lated approach; in the model the aircraft remains perfectly 
trimmed for a -3 deg glide slope approach. 

Fig. 5 presents solutions similar to the landing situa- 
tion of Fig. 4 except instead of sinusoidal gusts, the Boeing 
727 class airplane is experiencing horizontal wind shear of 
various magnitudes. There appears to be little significant 
difference in Ah' and Au* between headwind or tailwind shears. 
Both types of shears cause the airplane to stall somewhere 

along its landing flight path for shear values in the vicin- 

-2 

ity of + 4 m*sec 

Fig. 6 is the vertical wind shear equivalent of Fig. 5. 

For vertical downwinds which increase with time (positive 

values of w /t) the Ah' values are about the same magnitude 
U 

as the Ah' values for u ^/t (Fig. 5). However, vertical wind 

shears (+ w /t) do not produce values of Au ' as large as 
y 
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Fig. 5. RMS altitude and velocity deviations of a Boeing 
727 class airplane to horizontal wind shears during a land 
ing from 500 m altitude to the ground. 





UPWIND 


DOWNWIND 


Wg/t m sec"^ 

Fig. 6. RMS altitude and velocity deviations of a Boeing 727 
class airplane to vertical wind shears during a landing from 
500 m altitude to the ground. 




those produced by -h u g/t longitudinal wind shears. The large 

values of i^h ' produced by upwind shears are due to the fact 

that after a certain amount of time has passed, the vertical 

upwind is larger than the sink rate of the aircraft and the 

airplane never lands, but continues to climb. The /^h ' values 

for upwind shears have been computed for a time interval of 

133 seconds (the normal landing time for descent from 500 m 

with no wind) . All other ^^u ' and /i^h ' values shown in Figs. 4, 

5 and 6 have been computed by the following formulae; 

t^ 

Au' = [^ / u^dt]'^ (13) 

L-'o 

r^L % 

Ah' = (h-h^)2dt]' . (14) 

2. Aircraft Simulations Using Real Data 

The aircraft approach simulation model discussed in 
part III (1) has been used with actual measured winds. The 
two main sources for these winds have been 1) winds measured 
by research aircraft during flight, and 2) winds measured 
from a 444 meter-tall television tower. 

(a) NCAR Aircraft Data 

A prime source of data for winds near thunderstorms has 
been the National Center for Atmospheric Research (NCAR) 

Queen Air research aircraft. Such measurements were made, 
principally in the lowest 1 km of the atmosphere, on several 
days during the spring seasons of 1975 and 1977. These 
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aircraft are ideally suited for obtaining such data since 
they cruise at a speed close to that of the approach speed of 
a Boeing 727. The primary instrumentation system aboard the 
research aircraft is a gust probe wind measuring system cou- 
pled to an inertial navigation platform, providing 1.0 Hz 
sampling of the three-dimensional winds. Only the longitu- 
dinal and vertical wind components are used in this study. 

As mentioned, the data are 1.0 Hz samples of longitu- 
dinal and vertical winds. The conventions of positive values 
for tailwind and downwind have been used. Data from two days 
have been used in the simulations. These days are 13 June 
1975 (a day with tornadoes in Oklahoma) and 19 May 1977 (when 
the aircraft flew in the vicinity of a squall line and through 
the associated gust front) . About two hours of data are avail- 
able from each of the two flights. By initiating aircraft 
approach simulations at 15 sec steps into the data, a total 
of 480 simulations per flight have been produced. 

Figs. 7 and 8 give details of two sample simulations 
from the 13 June 1975 data. Similar simulations have been 
produced for other times that day and for 19 May 1977. The 
heavy straight lines show the -3 deg glide slope centerline 
and the so-called "missed approach" of -3 deg ^0.7 deg. 

The sequence of numbers "1" gives the actual aircraft posi- 
tion for each 2 sec from the initial time shown. The data 

panels give 1.0 Hz values of u , w , true airspeed (TAS) , 

^ y 
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Fig. 7. Fixed-stick model simulation using NCAR aircraft data 
for 13 June 1975 at 140230 CST. Data panels show u„, w„ input 
and TAS and 0 output, with time increasing to the right, with 
each division 10 sec. Larger panel shows aircraft position 
with respect to -3 deg +0.7 deg glide slope to nominal touch- 
down at lower right. Case demonstrates that seemingly turbu- 
lent input results in relatively good approach. 
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Fig. 8. Same as Fig. 7, for starting time of 142815 CST. 
Demonstrates non-turbulent but high shear input results in 
poor approach. 
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and pitch angle (0). The first case is interesting in that 
the approach is "quite good" since departures from nominal 
^TAS = 72.1 m* sec ^ (140 kt) and "on the beam"] are small, 
even though the input winds appear to be turbulent. The 
second case (Fig. 8) demonstrates a situation where the de- 
parture from the glide slope is extreme, while turbulent 
velocity deviations are relatively minor. 

To better quantify these "departures", variables Au' 
and Ak ' , as defined by Eqs . (13) and (14), were calculated 

for the 480 simulations for each day. Each quantity Au ' and 
Ah', can be considered as Approach Deterioration Parameter 
(ADP) under the assumption that significant velocity devia- 
tions from the nominal true airspeed and height deviations 
from the glide path represent deterioration of the approach. 
Table 2 gives the appropriate values for the example cases. 

Table 2. Approach Deterioration Parameters and shear values 
^o^ the two e xamples given in Figs. 7 and 8, 13 June 1975. 


Time 

Au* 

-1 

m- sec 

Ah' 

m 

(Ug/t)i 

-2 

m* sec 

^V^^30 

-2 

m* sec 

(Wg/t) 

-2 

m* sec 

( 30 
-2 

m- sec 

140230 

1.6 

14.8 

1.06 

0.078 

0.63 

0.040 

142815 

2.5 

193.7 

0.26 

0.176 

0.09 

0.017 


An attempt has been made to correlate shear values for 
u and w , to the values of ADP. Table 3 gives the least- 

g g 

squares linear regression correlation coefficients for four 
shear values. The subscripts indicate the niomber of seconds 
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for which the mean shear is calculated. Clearly for the Au' 
parameter, reasonably high correlations occur, while for Ah' 
no significant relationship appears, especially for the 13 
June 1975 data. The fact that Au' is highly correlated and 
Ah' is less highly correlated can be understood by examining 
the relationship between Au ' , Ah' and the shears as seen in 
Figs. 5 and 6. Au' is seen to plot as very nearly straight 
lines in Figs. 5 and 6 while Ab' does not, and is, in fact, 
a rather complex non-linear function. 


Table 3. Least-squares linear correlation coefficients be- 
tween indicated shears (independent variable) and two Approach 
Deterioration Parameters for 480 simulations (13 June 1975 
data). Values for 19 May 1977 data are shown in parentheses. 


Shear Values 


Au' 

(m-sec ^) 


Ah' 

(m) 



(0.63) 

(0.47) 


0.68 

0.09 


(0.50) 

(0.62) 

<^>30 

0.35 

0.30 


(0.47) 

(0.37) 

(w/t) 

0.69 

0.10 


(0.46) 

(0.43) 

< 30 

0.70 

0.01 


Variance spectra, as described in part II (3) have been 
found for each of the first 470 input data sets for the two 
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research days (13 June 1975 and 19 May 1977) . Each spectrum 
is calculated at 251 frequencies from 0.0 to 0.5 Hz. For 
frequencies at regular intervals, a least-squares fit is made 
correlating 470 pairs: the energy density at a given fre- 

quency and the value of Au ' for that simulation. A first 
order correlation coefficient is calculated for each fre- 
quency. Fig. 9 contains plots of the square of the correla- 
tion coefficient (positive correlation) against frequency for 
the two days. Two curves are shown in each plot, one for u 

9 

and one for w^. Correlations have been run for the special 
case of (.026 Hz) and are shown with the small x's. 

These results support the ideas presented in Fig. 4. 
That is, Au' is greater when a large amount of energy is pre- 
sent in the u^ field centered about the phugoid frequency. 

In other words, there is a strong linear relationship between 
Au ' and high energy density at ciUp^/ a relationship which 


breaks down at frequencies significantly different from u) 
However, a similar correlation peak is not seen in the w 


ph‘ 


field at u)pj^» nor, from Fig. 4, would it be particularly im- 
portant for aircraft approach quality were it present. Clear- 
ly then, we can conclude that longitudinal wind gusts provid- 
ing energy at the phugoid frequency may result in airspeed 
oscillations of a nature that would be difficult to control, 
and in fact may lead to stall, and otherwise disastrous re- 
sults . 


29 



(CORRELATION)" (CORRELATION) 




The Au ' values are plotted as functions of the peak value of 
the ARWS [discussed in part II (3)] in Fig. 10- Again both 
days are shown. The ARWS appears to be a better predictor 
for estimating approach quality. Close inspection of Fig. 10 
shows that ambiguities remain, since a precise Au' cannot be 
predicted from a given ARWS. Instead, only a critical thres- 
hold can be seen, beyond which large values of Au' might occur. 
An unanswered question is what value of ARWS is critical? 

This question will be addressed in part III (3) . 

(b) KTVY Tower Data 

Another source of data has been space-time adjusted 
fields of winds from the 444 meter instrumented KTVY-TV tower 
located in northern Oklahoma City, Oklahoma. Twenty thunder- 
storm influenced cases have been used to produce input winds 
for the simulation model. Approach wind conditions are ob- 
tained along approach paths in each direction through the two 
dimensional data fields. Also turbulence has been added by a 
numerical model of turbulence to create a second set of 40 
input winds for simulations. Thus a total of 80 approach 
wind inputs are available for use in the Boeing 727 simula- 
tion model. Tower data were processed by use of a computer 
program developed by Frost et al^. (1978) . 

Samples of these simulations are shown in Figs. 11 and 
12. These show extreme cases of aircraft performance. Fig. 

11 shows input winds with small scale turbulence but no large 
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Fig. 10. Scattergrams relating Au ' to ARWS 
function, expressed in dB, for (a) 13 June 
1975 and (b) 19 May 1977. JFK accident ex- 
ample shows a computed ARWS value of 65 m^. 
sec 2 *Hz 1 predicts a Au ' RMS velocity de- 
parture of 3.5 m*sec~l. 


32 








'5 


UG 

000020 




HORIZONTAL DISTANCE [KM] 


Fig. 12 . Same as Fig. 7 , for a sample simulation using tower 
data. Demonstrates wave-like input which results in a poor 
approach . 





fluctuations over sustained time periods. The aircraft per- 
formance is good with only slight airspeed changes. The air- 
craft does drift somewhat below the glide slope due to the 
downdraft (positive w^) . In Fig. 12, a large change in u^ 
over a number of seconds produces rather large fluctuations 
in airspeed and aircraft height. 

Fig. 13 is comparable to Fig. 9 except it is for the 80 
tower wind cases. Notice that the results are different. 

The correlations are smaller and the peak for the u^ spectra, 
near the phugoid frequency, is not seen. The smaller sample 
size (80 vs. 470) may be a factor in these differences. 

Fig. 14 includes scattergrams for Au' as a function of 
the peak value of ARWS . Three sets are included: 1) 40 cases 

of data with no turbulence added, 2) 40 cases with turbulence 
added, and 3) 80 cases consisting of the combination of both. 
The patterns are rather similar in each plot and agree with 
the larger sample shown in Fig. 10, and the conclusions reach- 
ed in the discussion of Fig. 10 apply here as well. 

3. An Example: JFK Crash, 1975 

The tragic crash of Eastern Airlines Flight 66 at JohnF. 
Kennedy Airport on 24 June 1975 has focused the attention of 
aeronautical engineers and meteorologists on the problem of 
landing aircraft in thunderstorm related environments. The 
National Transportation Safety Board in NTSB (1976a) has 
stated that the probable cause of Eastern Airlines ' Flight 66 
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Fig. 13. Same as Fig. 9, for eighty tower data simulations. 
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pig. 14. Plots same as Fig. 10, for (a) 
tower data cases, (b) tower data cases with 
turbulence added and (c) both (continued on 
next page) . 
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Fig. 14. Continued. 
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accident was "...an encounter with adverse winds associated 
with a very strong thunderstorm located astride the ILS local- 
izer course which resulted in a high descent rate into non- 
frangible approach light towers." 

We have examined the Eastern 66 crash using the Boeing 
727 simulation model. Although data are sparse for the winds 
encountered by the aircraft, the relatively simple descrip- 
tion of the winds found in NTSB (1976a) is used here. Fig. 15 


gives the u and w wind components taken directly from the 
9 9 

NTSB report (page 17) , and placed on a simulated -3 deg glide 


slope approach. Fig. 16(a) gives the simulated position of 

the aircraft relative to the ILS corridor, while 16(b) gives 

the simulated true airspeed. Three cases have been run: 

u with w = 0, w with u =0, and u w , all using the fixed- 
9 9 9 9 9 9 ^ 

stick mode. 


Several interesting observations can be made, related 
in particular to Fig. 4. When u^ is taken alone, input is 
an apparent one-half wave of period of approximately 19 sec, 
or a full wave of 38 sec; for such a wave the energy is con- 
centrated at 0.026 Hz, precisely the phugoid frequency! We 
give computed values of ADP as a function of frequency, so 
comparisons to Fig. 4 might be made; such a comparison for 

the u_ case is quite good. The case of the w fluctuation 
9 9 

must be considered in a little more detail. In Fig. 15, we 
see more of a pure shear of -10.8 m»sec~^ over about 12 sec. 
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Fig. 15. Longitudinal and vertical environment winds encountered 
by Eastern Flight 66, as derived from NTSB (1976a). Subjective 
smoothing was applied to the data; no turbulence fluctuations 
were added. A Ug headwind and a Wg updraft are positive. 
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Fig. 16. (a) Model output showing airplane position relative to 

glide path. Curves are shown for three input environmental winds 
longitudinal Ug, vertical Wg, and both Ug and Wg. (b) Model out- 
put showing airplane airspeed for three cases in (a) . Nominal 
airspeed is 72.1 m*sec"l (140 kt) . 




rather than a wave feature. To relate shear to harmonic fea- 


tures, a Fourier analysis of simple shear is performed, sug- 
gesting most of the energy density is represented in a sine 
wave with a 12 sec period. Consequently we can consider the 
appropriate frequency to be .083 Hz or uu/uUp^ = 3.21, as seen 
in Table 4. From Fig. 4, this frequency clearly represents a 
small deterioration of approach quality, as seen in the esti- 
mated Au ' and Ah'. Notice that the ADP ' s for = 3.21 in 

Fig. 4 are far below those calculated and given in Table 4. 
This difference may be due to the fact that we are comparing 
ADP estimates of two highly simplified models; however, the 
example serves to illustrate the point that the aircraft re- 
sponse to w is considerably less than its response to u . 

8 8 


Table 4. Approach Deterioration Parameters for JFK simula- 
tion (refer to Figs. 15 and 16 for profiles) . 


.... 


Au ' 

Ah' 

U) 


Case 

m* sec 

(m) 

(Hz) 



u 

g 

only 

3.5 

23.7 

.026 

1.00 

w 

g 

only 

0.9 

20.4 

.083 

3.21 


. w 

3.8 

22.6 




Returning to the full simulation, note that fluctuations 
in velocity from 72.1 m*sec ^ (140 kt) are large indeed, re- 
sulting in a wave-like departure from the nominal value. In 
terms of the Au' ADP, most of the deterioration occurs because 
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of the u factor. In fact, the calculated values of Au' and 
4h' agree quite well with the predicted values given in Fig. 4, 
for the appropriate value of We see that the largest 

values of ADP should be associated with u , and the most sen- 

g 

sitive indicator is the ^u' ADP. Referring again to Fig. 16, 
we see that the extreme airspeed oscillation associated with 
u bears this out. 

g 

The deleterious effect of w is present, as seen in the 

g 

glide slope position plot; however, we contend that the most 
critical feature is the extreme oscillation in airspeed, due 
to the aircraft's encounter with high energy at the phugoid 
frequency. 

Inspection of Fig. 4 indicates a large peaTc in the rela- 

_3 

tionship between w and Ah , near uu = 4.7 x 10 Hz represent- 

g 

ing a wave period of 214 sec. we assume that a real pilot 

easily would be able to overcome such an effect by elevator 

and thrust corrections, since the time-space scale length is 

61% longer than the 133 sec nominal approach period. Use of 

the fixed-stick mode precluded such a correction here. 

The utility of the ARWS function can be seen in the JFK 

simulation. For example, a Au' of 3.5 m*sec ^ requires the 

product of Ug energy density and transfer function be at 
2 -2 -1 

least 65 m -sec -Hz expressed xn dB ; as shown xn Fxg. 10. 
For the Eastern 66 simulation, this value of ARWS presiamably 
represents a critical threshold. 
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To obtain an estimate of the frequency of occurrence of 

conditions comparable to those experienced by Eastern 66, the 

u energy density at the phugoid frequency for the 940 near- 
9 

thunderstorm wind profiles was examined. The cumulative fre- 
quency for multiples of the energy density obtained from the 
wind profile used in the Eastern 66 simulation is shown in 
Table 5. It can be seen that for these near-thunder storm 
conditions, the energy density at the phugoid frequency ex- 
ceeds the Eastern 66 value for about 8% of the cases and is 
more than double that value for almost 2% of the cases. This 
demonstrates that flying conditions as bad as those experi- 
enced by Eastern 66 are not extremely rare. 

Table 5. Cumulative frequency distribution of energy density 
at phugoid frequency, as a function of multiples of Eastern 66 
computed energy density for Ug at (jupji* Tower data values are 
in parentheses. 

Multiples of Eastern 66 Energy Cumulative 

2 -2 -1 

Density of 95.24 m -sec -Hz Frequency 


0.25 

0.373 

(0.700) 

0.50 

0.203 

(0.363) 

0.75 

0.114 

(0.225) 

1.00 

0.078 

(0.200) 

1.25 

0.051 

(0.138) 

1.50 

0.029 

(0.125) 

1.75 

0.022 

(0.088) 

2.00 

0.017 

(0.050) 
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Included in Table 5 are the cumulative frequencies for 
the 80 tower data cases. It can be seen that 20% of these 
cases exceed the value of Eastern 66. The results are com- 
parable to the 940 aircraft cases but frequencies are some- 
what higher. 

In conclusion, we believe that the JFK accident is as- 
sociated with the airplane's encounter with a horizontal wind 
containing high energy at the airplane's critical phugoid fre 
quency, which caused a sudden extreme variation in the air- 
speed. Matters were made worse by the downdraft, but the air 
plane's response was most seriously affected by u^ variations 
We see this explanation as a subtle but important variation 
from the conclusion of Fujita and Caracena (1977) that a sud- 
den "downburst" drove Eastern 66 into the ground. Their ex- 
planation does not consider aircraft response stability as a 
vital parameter, and leaves the impression that the accident 
was caused entirely by meteorological events. Our explana- 
tion certainly includes meteorological events, including a 
possible downdraft as an initiating factor, but extends the 
concept to couple aircraft response factors to the weather. 
The horizontal wind u^, perhaps caused by thunderstorm out- 
flow, provides energy at the phugoid frequency, excites a 
severe horizontal velocity oscillation, leading to loss of 
control and a possible stall. Furthermore, we must not over- 
look the reality that an aircraft is likely to encounter 
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horizontal gustiness more frequently than vertical "downbursts" . 

At this point it is very interesting to compare Fujita's 
(1977) analysis of Eastern 66 to our model results seen in 
Fig. 16. Fig. 17 (from Fujita and Caracena, 1977) shows the 
flight track of Eastern 66, in a view similar to Fig. 16(a). 
There is a remarkable similarity between his analysis which 
requires the downburst explanation, and our u^ curve which 
includes no w^ component. True, then w^ is added, matters 
are made slightly worse, but the basis of our alternative 
explanation is made clear. We must conclude that a consid- 
eration of aircraft response theory, as well as environment 
wind energy, is necessary to fully appreciate the approach 
hazard . 
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Fig. 17. The path of Eastern 66 on 24 June 1975 in the vertical 
plane including the glide slope of runway 22-L at JFK (used with 
permission from Fujita and Caracena, 1977, Fig. 10). 



CHAPTER IV 


CONCLUSIONS 

The work reported herein represents an attempt to com- 
bine fundamental aircraft response to the thunderstorm envir 
onment wind shear, during an approach-to-landing. Not only 
have we developed several means of examining the nature of 
the wind shear, but we have mated it to the basic response 
characteristics of aircraft. This has resulted in a means 
to predict the possible deterioration of the approach, for 
given wind conditions and for a specific aircraft type. We 
have chosen only one aircraft type for our study, a medium- 
sized jet transport (Boeing 727 class) as a convenience, and 
not in any way to single out this airplane, but to represent 
a typical airliner. 

Our examination indicates the presence of wind shear 
containing energy at our aircraft's resonant or phugoid fre- 
quency sufficient to suggest serious deterioration of the 
approach on about 20% of approximately 1000 cases examined. 
This deterioration is observed as sudden oscillation in air- 
plane airspeed that can bring the aircraft near stall, or in 
height oscillations from the glide slope, possibly bringing 
about a premature impact with the ground short of the runway 

For our aircraft of study, the phugoid frequency is 
0.026 Hz, or for an approach speed of 72.1 m*sec ^ (140 kt) 
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a phugoid wavelength of 2.8 km. We are not surprised by our 
finding that approximately 20% of our cases had phugoid energy 
sufficient to suggest approach deterioration equal or worse 
than that encountered by E as tern_ Flight' 66. Such a scale 
length (near 3 km) would be a typical one in a thunderstorm, 
for vertical columns (updrafts and downdrafts) occur on such 
a horizontal length scale. Our study suggests that horizon- 
tal variations of the environmental wind along the aircraft ' s 
longitudinal axis of motion, are more critical than vertical 
variations of wind. While these vertical components can de- 
teriorate the approach, it appears that a firm handle on the 
horizontal component is most important. 

Consequently, we believe critical energy at the phugoid 
fr equency is at least not rare in the thunderstorm landing 
environment , and therefore the national airspace system should 
take steps to detect and warn approaching aircraft of such 
problems. It is important to realize that this conclusion 
about the nature of wind shear is fundamentally a time depen - 
dent problem. The existence of atmospheric waves, whose ef- 
fect on aircraft depends on approach speed as well, requires 
time variation considerations. Previous wind shear studies 
have not denied this important point, but have failed to 
recognize it as critical. 

Our work has limitations, to be sure. First of all, it 
is theoretical in that we simulate flight, rather than consider 
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a real aircraft. Furthermore, no pilot is considered. Some 
will suggest that our conclusions are severely limited by this 
exclusion. We cite our sensitivity studies of elevator input 
in chapter II and our Eastern 66 simulation in Chapter III as 
indicating that phugoidal response problems are still serious 
when pilot input is included. Finally, we have examined only 
one type of aircraft. While this of course limits the scope, 
our objective is not to identify positively whether a speci- 
fic airplane will in fact crash, but rather to indicate the 
fundamental nature of the wind shear - aircraft interaction. 
This interaction exists for all aircraft, but a given wind 
field may have less serious consequences for some aircraft 
(our Beech Baron at JFK, for example) . 

Where do we go from here? We believe a real-time ap- 
proach deterioration parameter detection system may be fea- 
sible. The system would work as follows: 

1. A single pulsed Doppler radar would collect radial 
velocities (from the ground) along intended ILS 
glide paths to airport runways. These radial com- 
ponents would be essentially identical to u^ com- 
ponents discussed in this work, after making a 
frozen turbulence hypothesis. Remember, our work 
has indicated the component should be a good 
"trouble" predictor. 

2. The u data would be collected at discrete distances 

g 
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along the glide path, then converted to the time 
series input form by knowing the aircraft approach 
speed. This u as a function of time would be fed 

g 

into the simulation model, which would be adjusted 
for a type of aircraft. 

3. The model output would predict an ARWS parameter 
which could be used by the ATC controller and/or 
pilot to decide whether severe conditions exist 
along the intended approach. 

The advantages of this system could be many. It would not 
require an actual and possible dangerous approach to be util- 
ized. A quantitative estimate of ADP would be available at 
all times. The ADP would be mated to aircraft type (although 
we do not expect the system to be used to let some aircraft 
to land, while forcing others to go around) . 

The authors are in the process of testing the feasi- 
bility of such a system, with the hopes that an operational 
system could be implemented in the early 1980 's. 
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